Stripes in the two-dimensional t-J model with infinite projected entangled-pair states 
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We simulate the t-J model in two dimensions by means of infinite projected entangled-pair states 
(iPEPS) generalized to arbitrary unit cells, finding results similar to those previously obtained by 
the density-matrix renormalization group (DMRG) for wide ladders. In particular, we show that 
states exhibiting stripes, i.e. a unidirectional modulation of hole-density and antiferromagnetic 
order with a 7r-phase shift between adjacent stripes, have a lower variational energy than uniform 
phases predicted by variational and fixed-node Monte Carlo simulations. For a fixed unit-cell size 
the energy per hole is minimized for a hole density pi ^ 0.5 per unit length of a stripe. The 
superconducting order parameter is maximal around 0.75 — 0.8. 

PACS numbers: 02.70.-c, 71.10.Fd, 71.10.Hf 



I. INTRODUCTION 

The simulation of strongly correlated electron systems 
in two dimensions remains one of the biggest challenges 
in condensed matter physics. The infamous negative sign 
problem^ prevents accurate quantum Monte Carlo sim- 
ulations of large fermionic systems at low temperature. 
An alternative route without a sign problem is provided 
by variational Monte Carlo (VMC). Typically, several 
ansatz wave functions, biased towards different orders, 
are optimized and the one with lowest energy is consid- 
ered as the best approximation to the ground state. A 
powerful way to improve upon a variational wave function 
is to use it as a guiding wave function for the fixed-node 
Monte Carlo (FNMC) method,!^ which yields the best 
variational wave function with the same nodal structure 
as the guiding wave function. 

In recent years, a new class of variational wave 
functions for two-dimensional systems have been pro- 
posed: the so-called tensor network states, including 
e.g. the multi-scale entanglement renormalization ansatz 
(MERA)"^ and projected entangled-pair states (PEPS),^ 
also called tensor product states (TPS).^ These tensor 
networks have recently been generalized to fermionic 
systems,^ and can be seen as generalizations of the 
density- matrix renormalization group (DMRG)^ method 
to two dimensions. In contrast to other variational wave 
functions, tensor network states are largely unbiased,^ 
with an accuracy which can be systematically controlled 
by the so-called bond dimension D (typically called m 
in DMRG). DMRG yields very accurate results for quasi 
one-dimensional ladder systems with cylindrical bound- 
ary conditions up to a width around 8 — 12, but becomes 
numerically inefficient for larger widths. The computa- 
tional cost of MERA and PEPS is polynomial in system 
size and in however, with such a high power in D 
that one may question if these methods are competitive 
to solve hard problems in condensed matter physics. 

In this paper we show that tensor networks indeed can 



compete with the best known variational methods, in 
particular with FNMC based on Gutzwiller projected 
ansatz wave functions.^ Specifically, we simulate the 
doped t-J model for J/t = 0.4 in the thermodynamic 
limit with infinite PEPS (iPEPS^^) and with DMRG for 
systems up to a width 10, and find significantly lower 
variational energies than obtained with FNMC. It was 
previ ously shown that iPEPS yields lower energies than 
VMC,™lbut the values were still higher than the ones 
from FNMC. Here we find that extending the ansatz to 
larger unit cells than the 2x2 cells previously used leads 
to a considerable improvement of the variational energy, 
which is an indication that the ground state may break 
translational invariance on a larger scale than 2 x 2.1^E^ 
By inspection of local order parameters we find that 
the ground state of the two-dimensional t-J model ex- 
hibits stripes, i.e. a unidirectional modulation of the 
hole density and the antiferromagnetic order, as previ- 
ously found with DMRG in cylinders up to width 8.1^ 
This is in contrast to the findings from VMC^^ and 
FNMC,^ which have suggested a uniform phase. 



II. METHOD 

Fermionic iPEPS has been introduced and explained 
in detail in Ref. |TT] and we here repeat only the basic 
ideas. A PEPS, illustrated in Fig.[l^), can be seen as an 
extension of a matrix-product state (MPS) in Fig. [TJd), 
the tensor network DMRG is based on, to two dimen- 
sions. Each blue-filled circle represents a tensor with a 
rank given by the number of legs attached to it, where 
the open leg corresponds to a physical index carrying the 
local Hilbert space of a lattice site, and the connecting 
lines are bond indices with a certain bond dimension m or 
which characterizes the number of variational param- 
eters in the ansatz. Tracing over all bond indices yields 
the coefficients of the state in the tensor product basis of 
the local Hilbert spaces of all sites. MPSs and PEPSs en- 



able an efficient representation of states obeying the area 
law of the entanglement entrop}'^ in one and two dimen- 
sions, respectively. An MPS can also be used to represent 
states in two dimensions, e.g. by using a snake structure 
as in Fig.[T]3), however, the required bond dimension m of 
the tensors grows exponentially with the system's width, 
whereas in a PEPS the required D is independent of sys- 
tem size (in the limit of large systems). In the present 
work we consider a PEPS of infinite size (iPEPS) made 
of a periodically repeated rectangular unit cell contain- 
ing Lx X Ly = N different tensors, A^^'^^ labeled by 
the coordinates relative to the unit cell. To obtain an 
approximate representation of the ground state we per- 
form an imaginary time evolution of an initial, randomly 
chosen, iPEPS. To efficiently compute observables an ap- 
proximate contraction scheme, discussed below, is used 
to evaluate the trace over all bond indices, where the 
accuracy can be controlled by another parameter called 
boundary dimension x- 

For the experts, we briefiy outline how to treat large 
unit cells, where we adopt the notation from Refs. ill|12[ 
First, for the im aginary time evolution we use the sim- 
ple updat^m^ on all bonds in the unit cell. Second, 
to contract the iPEPS we use the corner-transfer ma- 
trix (CTM) method^^ generalized to larger unit cells. 
The CTM method yields the so-called environment, con- 
sisting of corner tensors Ci,C2,C3,C4 and edge ten- 
sors Ti,T2,T3,T4, which account for the infinite sys- 
tem surrounding the tensors in the "bulk" of the system 
(Fig. [l]^)). We again assign coordinates [x^y] to each of 
these tensors to label the relative position in the unit 
cell, i.e. AN corner tensors and AN edge tensors in to- 
tal are separately stored. Initially, the corner and edge 
tensors at position [x^y] are constructed similarly as a 
reduced tensor a^^'^^, by multiplying the tensor ^[^'^1 to 
its conjugate and fusing the bond indices,'^ where we 
trace over the legs directed toward an open boundary, as 
illustrated in Fig. [T]p). The environment tensors are it- 
eratively built by four directional coarse-graining moves 
(left, right, top, bottom), similarly as proposed in Ref.l21l 
An entire CTM step consists of left moves, right 
moves, Ly top moves and finally Ly bottom moves. This 
sequence is repeated until convergence is reached. The 

is 



(a) PEPS 



(b) MPS 



renormalization procedure, adopted from Refs. I12|2C 
based on a 4 x 4 cell of tensors to compute an appropriate 
isometry. These isometrics are used to absorb a column 
(or row) of tensors into the corresponding boundary ten- 
sors, which effectively corresponds to the growth of the 
system by one column (or row) of sites. The left move 
is explained in Fig. [l] and one proceeds similarly for the 
other moves. 



III. SIMULATION RESULTS 

With the generalized iPEPS introduced in the last sec- 
tion we simulate the t-J model, given by the Hamiltonian 
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FIG. 1: (Color online) a) Diagrammatic representation of a 
PEPS with bond dimension D on a 4 x 4 lattice, b) An MPS 
with bond dimension m to represent a 4 x 4 system. c)-e) Dia- 
grams for the left move of the corner-transfer matrix method: 
(c) A plaquette of four reduced tensors a embedded in the 
environment. The coordinates [x^y] are to be understood rel- 
atively to the unit cell, i.e. one has to take the coordinate x 
(y) modulo Nx (modulo Ny). Cutting through the lines as 
marked in the figure yields the tensor Q^^'^^^^ in panel (d): a 
singular value decomposition of tensor q[^'^+^] is performed, 
where only the x largest singular values are kept. The re- 
sulting isometry [/[^'^+^] and its conjugate are used as an 
approximate resolution of the identity, I ^ f/t[a^,?/+i]f/[a^,y+i] 
For a fixed x = xo one computes all isometrics f/t^O'^] for 
all y G [1, A^y]. These isometrics are then used to obtain the 



,/[x-\-l,y] 



and 



renormalized corner tensors 
tensors 7^4^^^"'^'^^ for all 2/ G [1, A^^] and fixed x 



and edge 
xo, shown in 

e). This whole procedure is repeated Nx times for xq G [1, A^^] 
to complete an entire left move, (f) Initialization of a bound- 
ary tensor from a PEPS tensor and its conjugate, where cross- 
ings have been replaced by swap tensors (cf. Ref . [TT]) . 
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with a = {t, i} the spin index, hi = cJ^Qcr the elec- 
tron density and Si the spin 1/2 operator on site and 
Cia = Cia{l - c\^Ci^), whcrc we fix J/t = 0.4. 

Our main result for the infinite two-dimensional t-J 
model is that, in contrast to the uniform phase found by 
VMC and FNMC,^ we find a striped phase in the dop- 
ing regime 1/12 < 5 < 1/4, in agreement with DMRG 
calculations on wide ladders.^ Figure [2] shows that the 
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FIG. 2: (Color online) Energy per hole as a function of doping 
for J It — 0.4. For the D — % results the estimated error due 
to a finite x — 100 is 0.3%. For D = 10 the values are 
extrapolated in Grey arrows mark stripes with 0.5 holes 
per unit length. The DMRG cylinders labeled with a shift 
are wrapped periodically with the indicated shift to connect 
the transverse stripes into one long continuous spiral stripe, 
allowing arbitrary filling. 
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FIG. 4: (Color online) Examples of stripes running in the 
vertical direction obtained from various iPEPS simulations. 
Each panel shows one unit cell of the infinite lattice. The 
diameter of the dots scales with the local hole density with 
average values given by the (upper) red numbers. The arrows 
represent the local magnetic moment with average magnitude 
given by the (lower) black numbers. There is a vr-phase shift 
in the antiferromagnetic order between adjacent stripes. The 
width of a bond between two sites scales with the (singlet) 
pairing amplitude on the bond with a positive (green) or neg- 
ative (cyan) sign. A pattern with predominantly d-wave or- 
der is visible, with maximal pairing amplitude 0.01, 0.03, and 
0.003 in the three examples, respectively, a) Stripes of width 
4 in a unit cell 8x2 with hole density of ^ 0.5 holes per unit 
length per stripe at a doping ^ ~ 1/8. b) Stripes of width 6 
in a unit cell 12 x 2 with pi ^ 0.75, 3 ^ 1/8, where the pairing 
is maximal, c) Same as in b) but with pi — 1^ 5 where 
the pairing is suppressed. 



FIG. 3: (Color online) Energy per hole at doping 5 ^ 1/S 
for J /t — 0.4 as a function of truncation error in DMRG (left 
panel) and as a function of 1/D with iPEPS using a 8 x 2 unit 
cell (middle panel), with values obtained from extrapolation 
in the dimension x (right panel). The quadratic fit in the 
middle panel is a guide to the eye. 



variational energies obtained with iPEPS and DMRG are 
considerably lower than the ones from VMC and FNMC. 
We consider here the energy per hole (in units of t), 
Ehoie = {Eg — Eq)/5^ with Eg the energy per site and 
^0 = —0.467775 the value at zero doping taken from 
Ref. m The VMC and FNMC energies are seen to in- 
crease slightly with system size, and thus, we expect that 
in the thermodynamic limit the striped states are ener- 
getically favored. 

In Fig. |3]we present a convergence study of the iPEPS 
energies as a function of the dimensions D and Xi ^^^^ 
compare it to DMRG results for ladders for a doping 
6 ^ 1/8. While iPEPS is a variational ansatz the re- 
sulting energies may be non- variational, because of the 
approximate contraction of the iPEPS, which introduces 
an error depending on x- The convergence study in the 



right panel in Fig. [3] indicates that this error is smaller 
than the symbol size, i.e. the upper end of the symbol 
provides an upper bound of the ground state energy. The 
middle panel shows the values, extrapolated in as a 
function oil/D. The values have not yet converged, thus 
the ansatz can still be further improved by increasing D. 
Attempting a quadratic fit yields Ef^oie = —1.59(3) as 
D ^ oo^ which is similar to the extrapolated DMRG 
results, Ehoie = -1.61(2). 

Examples of stripes obtained with iPEPS are presented 
in Fig.[4j These stripes appear without biasing the initial 
iPEPS to stripe order (we typically start from several 
random initial states). We tested different unit cell sizes 
with Ly e {2,4} and G {2,4,5,6,7,8,10,12}, where 
each unit cell imposes a different periodicity on the wave 
function. Some of the Ly = 4 samples exhibit a slight 
modulation also in y direction, however, since the energy 
is similar to the Ly = 2 samples we focus on the latter. 

Figure |2] shows that, for a fixed unit cell size, the energy 
minimum is found near pi ~ 0.5 holes per unit length per 
stripe, e.g. at 5 = 1/8 for the 8x2 unit cell, which is in 
agreement with DMRG results (black stars and Ref. [16]) . 
However, the minimum for the 8x2 cell at 5 ~ 1/8 
is higher than the energy of the 10 x 2 or the 12 x 2 
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FIG. 5: (Color online) Mean pairing amplitude as a func- 
tion of linear hole density per unit lengtli of a stripe pi . The 
DMRG results are for one longitudinal stripe with a pairing 
field of 0.02 applied throughout, making the overall magni- 
tude somewhat arbitrary. 



sample at the same doping, which indicates that the re- 
pulsion between the stripes of width 4 in the 8x2 cell is 
too strong so that stripes with larger widths and larger 
pi are energetically favored, in contrast to the predic- 
tions by DMRG.^ For densities 0.75 < < 1 we ob- 
serve the tendency that mixed stripes (two stripes with 
different densities) yield a lower variational energy than 
two stripes of equal density, however, we do not observe 
signs of phase separation between stripes with pi = I and 
pi = 0.5 as found with DMRG^^. A future study with 
larger unit cells and larger D will shed further light in 
this issue. 

An open question is what happens at smaller doping 
than 1/12. The minimum at pi for different unit cells is 
seen to decrease with increasing unit cell length L^. It 
is conceivable that this trend continues, so that with de- 
creasing doping we obtain stripes which are increasingly 
more widely spaced, with a distance between the stripes 
varying as d ~ 1/{2S)}^ However, from the present data 
we cannot rule out other phases, such as a uniform phase 
at small doping S < 1/12, or phase separation between 
an undoped and a doped region. 



Finally, we study the pair field A 



l/V2{ciCj 



between nearest-neighbor sites (z, j) which is modulated 
along the x-direction and predominantly forms a d-wave 
pattern as shown in Fig. |4] The d-wave order between 
neighboring stripes has the same phase, however, we have 
also observed states with a similar energy where the d- 
wave order exhibits a 7r-phase shift between neighbor- 
ing stripes. Thus, it seems that this 7r-phase shift has 



only little influence on the energy, and from the present 
data we cannot determine which state is preferred. The 
mean pairing amplitude (averaged over the unit cell) in 
Fig. [5] exhibits a maximum around pi = 0.75 (p^ = 0.8 
in DMRG). For pi = 1 the stripes are insulating with 
a vanishing pairing amplitude. For pi = 0.5 the pair- 
ing amplitude for D = 8 is finite, but decreases quickly 
with increasing D, and is possibly insulating for larger 
too. In between these two insulating states, excess 
holes (or electrons) form pairs, leading to an increase of 
the pairing. 

IV. CONCLUSION 



In this paper we presented iPEPS simulation results 
for the t-J model in the thermodynamic limit, where 
we obtained striped states, in qualitative agreement with 
DMRG calculations on wide ladders. The iPEPS vari- 
ational energies are considerably lower than the ones 
from the uniform states obtained with VMC and FNMC 
based on Gutzwiller projected ansatz wave functions, 
and compatible with DMRG results for finite systems. 
This demonstrates that iPEPS is a competitive varia- 
tional method for strongly correlated electron systems.!^ 
The differences between iPEPS and DMRG only involve 
quantitative details, such as the precise stripe linear fill- 
ing and spacing as a function of overall filling. Using 
larger unit cells in iPEPS than the usual 2 x 2 is es- 
sential to study ground states which break translational 
invariance on a larger scale than only 2 lattice sites. 

The iPEPS wave functions for the accessible values of 
D include unbiased quantum fluctuations over short dis- 
tance and high energy scales. Over longer length and 
lower energy scales, they revert to a more mean- field-like 
description. Thus, with the currently accessible I^'s^^, 
we cannot expect them to resolve between, e.g., static 
versus fluctuating stripes. On the other hand, our re- 
sults provide significant evidence that approaches that 
ignore stripes do not give reliable descriptions of the t-J 
model. 
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